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f\ ' yield some first insights into the strong cosmic censorship issue and the conjectures 
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the current investigations suggest that strong cosmic censorship holds in this 

class. We further identify open issues in our current approach and point to future 

research projects. 
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1. Introduction 

Motivated by the desire to understand the history and fate of our universe, studies of 
cosmological solutions of Einstein's field equations have a long tradition. Observations 
indicate that there was a big bang in the distant past, and indeed, the simplest 
cosmological models, namely the Friedmann solutions for reasonable matter fields, 
predict such a collapse of spacetime. Now, the question arises whether such curvature 
singularities occur generically in solutions of Einstein's field equations or only under 
the strong symmetry assumptions underlying the Friedmann models. This problem 
was solved, namely due to the Hawking-Penrose singularity theorems |2B] we can 
expect incompleteness of causal geodesies in a wide class of situations. However, the 
information about the origin of the singular behavior provided by these theorems 
is very limited. One would expect that curvature singularities develop, as in the 
Friedmann case, which 'stop' the geodesies. But one was able to find solutions, as 
the A-Taub-NUT spacetimes discussed below, where the curvature stays bounded 
but the spacetime becomes acausal in a well-defined sense. Let us restrict in all 
our investigations to vacuum with positive cosmological constant A, in which case 
Einstein's field equations (EFE) in geometric units c — 1 and G — 1/{8tt) read 

Rfii/ — Xg^u, (1) 

and to four spacetime dimensions. Here, g^^^ is a Lorentzian metric of signature 
(—,+,+,+); Rfiu represents the Ricci tensor of g^^. One particularly important 
prototype solution of these equations is the de-Sitter spacetime [26] . The positive 
cosmological constant is motivated by current observations [43j |44j . Under these (but 
actually under much wider) assumptions, the Cauchy problem of EFE is well-posed 
and leads to the notion of the maximal globally hyperbolic development (MGHD) of 
a given Cauchy data set. See |2] and references therein for a review of the details. 
Now, for Cauchy data corresponding to a A-Taub-NUT spacetime for instance, the 
corresponding MGHD is extendible in several incquivalent ways as smooth solutions 
of EFE. In this sense, Einstein's theory has only limited predictive power. The hope, 
expressed by strong cosmic censorship (SCC) essentially first formulated by Penrose 
[35l[36], is that, although examples like the A-Taub-NUT spacetimes exist, the MGHD 
corresponding to a generic Cauchy data set is inextendible. This would rescue the 
predictive power of Einstein's theory. So far, there is neither a clear idea which 
topology on the space of initial data is supposed to give rise to the notion of 'genericy' 
here, nor which class of extensions to choose in general. Although this issue has only 
been solved successfully so far in special situations, some of which are commented on 
below, interesting new insights have been obtained as summarized in [21 I39j . 

Besides strong cosmic censorship there are other open problems of interest for 
cosmological spacetimes. Belinskii, Khalatnikov, and Lifschitz [301 El E] conjectured, 
what we nowadays call BKL-conjecture, a description of the properties of solutions 
with curvature singularities. We give only few further details later, and just remark 
that it is still a matter of active research whether their conjecture holds generically 
[21 [39] . Another outstanding problem is to find a characterization of the development 
of generic cosmological solutions in the expanding time direction, and the notion of 
cosmic no-hair essentially introduced in [27] is one promising scenario. 

In our work, we focus on the class of Gowdy spacetimes with spatial S^-topology 
which yields one of the main technical complications. In that class, the issue of SCC 
and the BKL-conjecture is still unsolved; however, there are strong expectations due 
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to results in the Gowdy case with spatial T'^-topology. We restrict to a special class of 
solutions which obey a generalized notion of cosmic no-hair by construction. Further, 
this class is close to A-Taub-NUT in the sense that it incorporates a "perturbation" 
parameter e such that e = corresponds to a A-Taub-NUT spacetime. Our solutions 
are constructed using Friedrich's notion of the J''*'-Cauchy problem discussed below 
by means of a numerical code introduced in [9] . The purpose of our investigations is 
to check whether SCO is true in this class and how the dynamics gives rise to it when 
e approaches zero. We are further interested in the BKL-conjecture. We see these 
investigations in this special class as first steps in order to shed further light on the 
general SCC and BKL-conjectures. 

The paper is organized as follows. First we give the relevant background material 
in section [21 including some remarks about the conformal field equations and future 
asymptotically de-Sitter spacetimes, the definition of the A-Taub-NUT family, a short 
introduction into the class of Gowdy spacetimes with spatial S'^-topology, our family 
of initial data and a few comments about our numerical method. In section [31 we 
write our main result, support it by means of numerical computations and heuristic 
arguments and discuss the current limitations of the approach. Finally we conclude. 

2. Background material 

2.1. Notations and conventions 

In this paper, we will always assume Einstein's summation convention when we write 
tensorial expressions. Any tensor can be represented by the abstract symbol, say, T 
or by abstract index notation, e.g. T^^^ , dependent on the context. Note however, that 
when we write such an indexed object, we can mean either the abstract tensor or a 
special component with respect to some basis. Hopefully this will be always clear from 
the context. Our convention is that Greek indices ^,iy, . . . = 0, . . . , 3 refer to some 
choice of local spacetime coordinates, whereas Latin indices i, j, . . . — 0, . . .3 represent 
indices with respect to some orthonormal frame field. When we have chosen a time 
coordinate x^, then Greek indices a, (3, . . . = 1,2,3 represent spatial coordinates, and 
for a choice of frame {ci} with timelike vector field eo, Latin indices a,b, . . . = 1,2,3 
refer to spatial frame indices. Writing {ci} just means the collection of tangent vector 
fields {eg, . . . , 63}. If a 2-indexed object is written in brackets, for instance (T^^), we 
mean the matrix of its components, where the first index labels the lines of the matrix 
and the second one the columns. 

2. 2. Future asymptotically de-Sitter spacetimes and the conformal field equations 

Let a smooth Lorentzian manifold {M^g^y) be given; not necessarily a solution 
of EFE. According to Penrose [34l [36], it is said to have a smooth conformal 
compactification if there exists a smooth Lorentz manifold-with-boundary (M, g) with 
boundary J := dM and a smooth function fi : Af ^ R such that 



(i) there is a diffeomorphism $ : M ^ M\J such that g = $» I J7 ^^Lrv ^ 

(ii) we have fJ > in the interior of M, and, $1 = and dO, 7^ on J . 

For simplicity we will always identify M with $(M) C M . Under these assumptions, 
J7 is called conformal boundary and can be interpreted as representing "infinity" of 
the physical spacetime (M,g^y). In the following, quantities which arc related to 
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the physical metric g^j^^ are marked with a tilde, while those related to the conformal 
metric 5^,^ are not marked. For a recent discussion and several examples see [T^. The 
main prototype solution of EFE with A > with a smooth conformal compactification 
is the de-Sitter solution. Further examples are provided by the A-Taub-NUT family 
discussed below. Now, in general if {M,g^^) is a solution of (H]) in vacuum with 
A > 0, then the corresponding conformal boundary J^ must be spacelike. In this 
case, we define J^ as the intersection of J with the chronological future of M in M 
and analogously J~ . We assume in all what follows that J is the disjoint union of 
J'^ and J~ , each of which is allowed to be empty so far. If J^ is non-empty, we 
call {M,gf^„) future asymptotically de-Sitter (FAdS); in the same way one can define 
past asymptotically de-Sitter solutions. According to [3], a FAdS spacetime with 
compact J^~^ is globally hyperbolic if and only if it is future asymptotically simple, 
i.e. all future directed inextendible null geodesies starting in M have future endpoints 
on J'~^. Then, all Cauchy surfaces of {M,g^^) are actually homeomorphic to J^ . 
From the physical point of view it is further important that all FAdS solutions obey 
a generalized notion of the cosmic no-hair picture [7] , and hence provide a model of 
the current observational predictions for the future of our universe. 

In order to compute solutions of the EFE in terms of conformal variables 
directly, Friedrich succeeded in finding a formulation of EFE called conformal field 
equations (CFE), which is both conformally invariant and regular even on J. By 
conformal invariance we mean that in addition to the coordinate gauge freedom of 
Einstein's theory, the CFE obey the conformal gauge freedom; i.e. the possibility 
to conformally rescale the conformal metric with a smooth positive function leaving 
the physical metric invariant. By 'regular' we mean that these equations provide 
evolution equations which are hyperbolic even in J". A review of his results and 
further references are given in [19]. His formalism is restricted to four spacetime 
dimensions, the relevant case for our investigations. Anderson [1] was able to find a 
different formulation applicable in arbitrary even spacetime dimensions. We, however, 
restrict here to Friedrich's CFE which we write in a special form below. Under the 
assumptions above, in particular for A > 0, the CFE allow one to formulate what we 
call JT """-Cauchy problem [TB]; Anderson's generalization of this concept was proposed 
in [1]. The idea is to prescribe data for the CFE on the hypersurface J'~^ including 
its manifold structure subject to certain constraints implied by the CFE. These data 
can be integrated into the past by means of evolution equations implied by the CFE. 
The J'"'"-Cauchy problem is well-posed and the unique FAdS solution corresponding 
to a given choice of smooth data on J"^ is smooth as long as it can be extended into 
the past. It is remarkable that this Cauchy problem allows us to control the future 
asymptotics of the solutions explicitly by deciding on the data on J'^ . Concerning 
the past behavior, the concepts of MGHDs and SCC with respect to data on JT""*" 
are well defined, however, there is only limited a priori control in the same way as 
it is the case for the Cauchy problem of EFE. In this paper, we will give no details 
on how the constraints implied on J^^ look and say only briefly what the relevant 
data components are, since we do not want to introduce all the geometric concepts 
now. Instead we refer to [El [7] where the details have been carried out. In order to 
choose a complete set of smooth data on JT""*", besides certain quantities related to the 
gauge, it turns out to be sufficient to choose the differentiable manifold structure of 
J'~^, which we assume to be the standard one on S"^, the induced conformal 3-metric 
hafj and the electric part of the rescaled conformal Weyl tensor Eafs [19l [20] . From 
those components the complete data set is fixed in terms of algebraic or differential 
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relations. The choice of these components is unconstrained except for E^fs which has 
to satisfy D°'Eap = 0. Here D is the Levi-Civita covariant derivative operator of /i„^, 
and the index in the constraint has been raised by means of hap. Befow, we mention 
a special class of solutions of this constraint. 

For our purposes, we decided to use the so-called general conformal field equations 
which are the CFE in conformal Gauf5 gauge [HI [19] . This fixes the gauge freedom for 
the coordinates, the conformal factor and a frame which is orthonormal with respect to 
the conformal metric. In our applications, we specialize the gauge even further to what 
we call Levi-Civita conformal Gaufi gauge [7]. One assumption for that is that the 
conformal geodesies representing the timelines leave J^'^ orthogonally. In particular 
this implies that they correspond to physical geodesies up to parametrization |22j . 
Hence this gauge can be considered as a (standard) Gauf5 gauge with respect to the 
conformal (and also the physical) metric up to parametrization. We can expect that 
this simple choice of gauge can lead to serious problems in the evolution, and in fact, 
we will discuss problems related to the gauge later. In any case, assuming without 
loss of generality A = 3 and having fixed the residual gauge initial data as described 
in [7], a hyperbolic reduction of the general conformal field equations in Levi-Civita 
conformal Gaufi gauge is 



dtel = -xy,\ 




(2a) 


dtXab = -XaXcb - ^Eab + -^afcj 




(26) 


dtr,\ = -x/r/, + nB,de\', 




(2 c) 


dtLab = -dtfl Eab ~ Xa^Lcb, 




(2d) 


dtEf, - De.i?a(/e''^e) = -^Xc'Efe + ^X^^Ejy - 


- Xc'E,'9ef, 


(2e) 


dtBf, + i?e.i?a(/e"'e) = ''^XcBfe + 3X(/S/)c " 


- Xc'B.'gef, 


m 


m) -^t(2-i), 




(25) 


for the unknowns 






u^{e,\xab,T^\,Lab,Ef,,Bf,). 







Among the unknowns are the spatial components ej' - with respect to some reference 
basis specified below - of a smooth frame field {e^} which is orthonormal with respect 
to the conformal metric and where, in this special gauge, eo = dt is the past directed 
timelike unit normal of the t — coTisi-hypersurfaces, the spatial frame components 
of the second fundamental form Xab of the t = const-hypersurfaces with respect to 
eo, the spatial connection coefficients r^**^, given by ^J'c^b = Ve^Cc — Xac^o where V 
is the Levi-Civita covariant derivative operator of the conformal metric, the spatial 
frame components of the Schouton tensor Lab, which is related to the Ricci tensor of 
the conformal metric by 

L^„ — Rfj.iy/2 — gfj_^g'"^R per/ 12, 

and the spatial frame components of the electric and magnetic parts of the rescaled 
conformal Weyl tensor Eab and Bab defined with respect to eo. In this special 
conformal Gaufi gauge, the timelike frame field eo is hypersurface orthogonal, i.e. 
Xab is symmetric. In order to avoid confusions, we point out that the conformal factor 
fl is part of the unknowns in Friedrich's formulation of the CFE. However, in the 
special case of vacuum with A > 0, it is possible to integrate its evolution equation 
in any conformal Gaufi gauge explicitly [H] so that it takes the special form ( [2g| in 
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our gauge. Hence, if the solution develops a smooth J^^ then it must correspond 
to the t = 2-hypersurface. The initial hypersurface J''^ corresponds to t — 0, and 
t increases into the past. We note further that, since Eat and Bab are tracefree, we 
identify E33 — —En — E22', the same for the magnetic part. The evolution equations 
of Eab and Bab are derived from the Bianchi system [19]. In our gauge, the constraint 
equations implied by the Bianchi system take the form 

D,^E'^ - r,-\BaaXb'' = 0, D,^B\^ + v'^KEdaXb" = 0, (3) 

where "q is the volume form of the spatial conformal metric with 77123 — 1, and 
indices are shifted by means of the conformal metric. The other constraints of the 
system above are equally important but are ignored for the presentation here. Further 
discussions of that evolution system and the quantities involved can be found in the 
references above. 

As the reference basis tangential to the t = consi-hypersurfaces assumed to be S'^ 
topologically, we choose frame vector fields {Ya} which are SU(2)-left invariant and 
satisfy the commutator relations 

3 

[YaM='^Y.'labcYc^ 
c=\ 

where Tiabc ^^ ^^^ totally antisymmetric symbol with 77^23 ~ 1- These requirements 
fix those vector fields uniquely up to the specific choice of realization of the SU(2)- 
action on S^ and up to certain transformations on the three-dimensional real vector 
space of SU(2)-left invariant vector fields. We do not list our choices involved here but 
rather write the explicit expressions belo"v4B- Henceforth we write Ca = e^^'YI,, and in 
the Bianchi system ((2e|) and ^^ and corresponding constraints (l3|), the orthonormal 
frame fields are interpreted as linear differential operators whose coefficients are given 
by Ca and by the components of \Ya\ in some coordinate system which we fix in a 
moment. 

We finally note how the Kretschmann scalar of the physical metric can be 
computed in terms of the variables above if the physical metric is a solution of EFE 
in vacuum with A > 0, 

~K = K^'^'h.u = 24y + m"" J2 d^afcl' - \Bab\^) ■ (4) 

a,b=l 

2.3. X-Taub-NUT solutions 

Consider the function F : R ^ R defined by 

V{t) = DoT^ + 2{iDQ - 2)t2 + Cor + 4 - 3i?o 
where [Dq^Cq) G M'^° x R are arbitrary parameters. For given (Do,Co), let / be 
a maximal interval on which y > 0. Corresponding to the choice of cosmological 
constant A > 0, parameters (DqjC'o) and interval /, we define the A-Taub spacetime 
[TT] as the Lorentz manifold (M, g) with M = / x §^ and 

»4"» 



r i + t2 

V V{t) 


dT^ 


, V{t) 


cr^(g)cr^ 


+ (l + r 


^)(- 


^ «)Cr^ + 


J^^a^)" 



I To be more precise we choose {Ya} as the reference frame on the spacetime as follows. Let {Ya} 
be given as before on one t = const Cauchy surface. Then we extend the frame to the spacetime, 
keeping the same symbol {Ya}, by means of the conditions [dt, Ya\ = 0. 
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which can be shown to be a solution of Rfj,i, — \g^v Here, {cr°} is the dual of the 
frame {Ya}- Note that for a given choice of (DqjCo) there can exist several distinct 
choices of the interval / and hence several distinct A-Taub spacetimes. 

The A-Taub spacetimes are globally hyper- 
bolic and maximal with this property. The de- 
Sitter solution is recovered by choosing the pa- 
rameters (Do, Co) = (A/3,0), in which case / = 
R. Indeed, due to the nonlinear stability of the 
de-Sitter spacetime [I71[T], there is an open neigh- 
borhood U of the de-Sitter point in the param- 
eter space R^° x R in which all corresponding 
A-Taub spacetimes are causal geodesically com- 
plete and asymptotically de-Sitter in the future 
and in the past. In general, for any given choice 
of (Z?o, Co), one can find tq e MU{— 00} such that 
one pick / =]to,oo[, and we will assume such I 
in the following. The corresponding A-Taub so- 
lution is FAdS, as depicted in figure [TJ note that 
in our convention d^ is future directed. Now, if 

V has zeros on R, i.e. if tq ^ —00, we will say 
that the corresponding A-Taub spacetime is sin- 
gular. Indeed, the metric is clearly singular when 

V becomes zero, however, curvature invariants 
are bounded. One can show analogously to the 
case treated in [T3] that singular A-Taub space- 
times are always, independent of the type of zero 

of V dX Tq, extendible and causal geodesically incomplete. In particular, this provides 
counter examples to theorem 1.3 in [T]. Any extension of a singular A-Taub solution 
is non-globally hyperbolic and will be referred to as A-Taub-NUT. We will not write 
extensions here but only note the following. The function V can only have a single 
or double real zero at r = tq. In both cases the r = ro-hypersurface corresponds to 
a smooth Cauchy horizon diffeomorphic to S'^ ruled by closed null curves. In the first 
case we will call the Cauchy horizon "non-degenerate" , in the second case "degener- 
ate" . This corresponds to the terminology of [22 which can be seen by constructing 
extensions in "Gaussian null coordinates" . figure [1] depicts a conformal diagram of a 
singular A-Taub solution. In that figure, circles correspond to r = consi-hypersurfaces, 
the bold circle represents J^ and the dashed circle the Cauchy horizon. We will be 
particularly interested in the non-degenerate case here. Note that we will often abuse 
terminology and speak, instead of A-Taub, often of A-Taub-NUT. 

Singular A-Taub solutions seem to be counter examples to SCC. Within the 
Bianchi IX class for A = 0, the analogous Taub-NUT spaces [ISl [55] for A = 
have been shown to be non-generic in a well-defined sense [40]. For A > one 
believes that the same holds, but no such rigorous result has been found yet. A more 
general extensively studied symmetry class is Gowdy symmetry which we discuss in 
the following subsection. In general, further relevant rigorous results concerning SCC 
for solutions with spatial S'^ -topology are [321 HIl E] • 



Figure 1. Causal 
structure of a singu- 
lar A-Taub solution. 
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2.4- Coordinates and Gowdy spacetimes 

Gowdy symmetry [25 is defined by the presence of a smooth effective isometric action 
of the group U(l) x U(l) with spacelike orbits. As summarized in [T3], if we restrict to 
smooth orientable connected compact spatial 3-manifolds, the only admitted spatial 
topologies are T^, S^ x §^ and S"^ (and lens spaces which are always included in the 
following discussion without further comment). Any two of such group actions on one 
of these topologies are equivalent in the sense of [ISj . 

For the latest state of knowledge in the class of T'^-Gowdy spacetimes, we refer to 
[411 142] and references therein; we rather focus on the S'^-case here. A representation 
of the action of the Gowdy group on S'^ is given by 

y y 

xi = COS — cosAi, X2 = COS — sinAi, 

X X ^ ' 

x^ — sin — cos A2 , X4 — sin — sin A2 , 

where we assume the standard embedding of §'' into M.*. The parameter x G [0;'"'] 
labels the group orbits while Ai and A2, both defined on Mmod27r, can be considered 
as the group parameters. The action degenerates exactly when y = or x = tt, 
i.e. the group orbits become one dimensional. These submanifolds of §^-topology are 
referred to as symmetry axes in the following. By abusing the terminology slightly, we 
often also speak of the south pole when we mean the symmetry axis given by x = tt; 
analogously, by north pole we mean the axis given by x — 0- Equation ([5]) can also 
be considered as the defining equation for a coordinate patch on S'^; these coordinates 
(x, Ai,A2) are called Euler angle parametrization. We will however rather use the 
coordinates (X;Pi:P2) determined by 

Ai = (pi +p2)/2, A2 = (pi - p2)/2. 

The Euler angle parametrization covers smoothly the dense subset of S"^ C §'^ which is 
S^ minus those points given by x = or x = ■"■ where the coordinates become singular. 
Our choice of representation of the frame {Ya} yields the following expressions 

Yi = 2 sin pid^ + 2 cos pi (cot x dp^ — esc x^p^) ^ (6a) 

I2 = 2 cos pid^ — 2 sin pi (cot x dp^ — esc xdp^) , (66) 

^3-29,,. (6c) 

Now, on a Gowdy-invariant spacetime with spatial S'^-topology let us choose 
spacetime coordinates (i,x, Ai,A2) such that the t — const-hypersurfaces are Cauchy 
surfaces, and such that on a i = io-hypersurface Eq the Gowdy orbits are x = const- 
surfaces parametrized canonically as above by Ai and A2. On Eo, a Killing algebra 
basis is {dxi^dx^}, but also {1^3,^3} where {Za} are SU(2)-right invariant vector 
fields defined in the same way as {Ya} before; in particular Z3 = 2dp^. Let us now 
extend the fields Ys and Z3 from So to the spacetime by means of the conditions 
[dt, I3] = \dt, Z3] = 0. Then it is clear directly that {Y3, Z^} is a basis of the Killing 
algebra at every instance of time if and only if this is the case for {dx-^jdx^}. It is 
obvious that the latter does not need to be fulfilled for all foliations; however, it is 
true in our gauge [J and also in areal gauge defined below. 

The issues of SCC and the BKL-conjecture in the S^-Gowdy class have so far 
only been solved in the polarized vacuum case with A = [221 El]- In particular, there 
is a nice characterization of those non-generic Gowdy-invariant data whose MGHD 
is extendible through a Cauchy horizon. For the non-polarized case, both issues are 
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open. From the results obtained in the T'^-case, one expects that spikes are present 
generically at the singularity, but in particular the behavior at the symmetry axes is 
not completely understood [35] . The purpose of our investigations is to shed further 
light on sec and the issue of spikes for §^-Gowdy perturbations of the A-Taub-NUT 
family. Similar numerical investigations in the §^ x §^-case have been carried out in 
[24]. 

Let us introduce the notion of the orbit area density. This quantity is defined as 
y^'iet{g{d\^,dx^)) (a ~ 1,2). In the S'^-case one can show that it is proportional to 
sin X and we will speak of the rescaled orbit area density when this factor has been 
divided out. For A = in vacuum, Chrusciel [T3j succeeded in showing that an open 
subset of the MGHD corresponding to generic Gowdy-invariant Cauchy data can be 
covered by an areal foliation, i.e. a foliation of spacelike hypersurfaces given by the 
level sets of the rescaled orbit area density, and that the orbits reach zero area either 
at the boundary of the MGHD or in the interior. In the polarized case, the result 
is stronger, namely, that the boundary of the MGHD is attained exactly where the 
orbit area vanishes in the generic case. There are no such results for A > 0. In 
particular, it cannot be expected that solutions which collapse and later re-expand 
due to the repelling force of the cosmological constant can be covered completely by 
areal foliations. 

Other aspects of Gowdy symmetry and the more general U(l)-symmetry which 
are needed to develop our numerical method are summarized in [9]. In particular, 
it is explained there how our formulation of the equations can be reduced to (1 + 1) 
dimensions in the Gowdy case in an indirect manner, since it is not possible directly. 
We should point out that this reduction allows us to restrict to the coordinate locations 
pi = and an arbitrary (fixed) value of p2 ■ Hence effectively the unknowns depend 
only on (i, x) and our problem reduces to 1+1 in the Gowdy case. All our computations 
here take place at pi = 0; the solution at other points can be computed from the 
formulae in [3] straight forwardly. 

2. 5. Class of initial data 

As initial data on »/+ for the CFE in Levi-Civita conformal Gaufi gauge, we use the 
data derived for U(l)- and Gowdy symmetry in [Tj, see also ^9 . Here we restrict to 
Gowdy symmetry. We make no claims about the genericy of these data even in the 
class of Gowdy data. Under the conventions above, the first part of the data take 
following form 

(e,") -diag(l,l,a), (7a) 

(Xafc) =diag(-l,-l,-l), (76) 

r/2 =0, T,\ = ~l/a, T^\^0, T^\^l/a, (7c) 

^2% =0, r,\^l/a-2a, Fgig^O, T,\^0, [Id) 

(Lab) = diag((5 - 3/a2)/2, (5 - 3/a^)/2, (-3 + 5/a^)/2^ , (7e) 

(Bab) = diag(~4(l - a^)/a^ -4(1 - a^)/a^, 8(1 - a^)^^) . {7f) 

The induced conformal 3-metric of »/+ is a Berger sphere with anisotropy parameter 
a > and an adapted orthonormal frame which is expressed by (|7a|) . Now, the only 
inhomogeneous part of the initial data is given by the components of Eab- In order 
to construct solutions of the constraint D"-Eab = on the background defined by 
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the expressions above, we expand the components of Eat in terms of spin-weighted 
spherical harmonics on S^, which under the assumption of U(l)-symmetry are related 
to spherical harmonics Wnp on §^ . For our definition of the functions Wnp consult [7] . 
We restrict to solutions with n = and n = 2. Then the expressions are 

(Ef)+ew2o — v^aeRew2i 

Eo + ew2o -^/2aeImw2l 

— V2aeRew2i —\/2asIiaw2i — 2(£'o + £^20) 
With respect to the Euler angle coordinates at pi = we have 
W20 = cosx, W21 = sinx/%/2. 

The parameter Eq is associated with the part of the data given by n = 0, while e 
parametrizes the inhomogeneous n — 2 part. 

Note that the solution corresponding to a general choice of parameters in this 
family is of unpolarized Gowdy type. This has already been claimed in [7J; however, 
there the argument turned out to be incomplete. In general, we call two initial data 
sets equivalent if the corresponding MGHDs are isometric. In particular, two initial 
data sets on J^^ are equivalent if there is a diffeomorphism from ^7+ to itself such 
that the puUback of all data fields equals the corresponding fields in the other data 
set, up to gauge transformations. In our special class of initial data now we find that 
two initial data sets, one given by £1 and the other by £2 = — £1, both with the same 
a > and Eq G K, are equivalent by means of the diffeomorphism from §^ to itself 
represented by {xi,X2,X3,X4) t—f {x3,X4,Xi^X2). In terms of the Euler coordinates, 
this diffeomorphism takes the form (x, /0i,p2) 1-^ (tt — XjPIj ~P2)- Hence without loss 
of generality we can assume that £ > 0. For £ = 0, the corresponding MGHD of J'^ is 
a A-Taub spacetime whose parameters {Dq, Cq) can be identified with our initial data 
parameters by 

a^l/^n, En ^ Co^/d'o/2, 

for A = 3. In Levi-Civita conformal Gaufi gauge, all r-level sets are also i-level sets. 
Note however that our convention is that dr is future directed while dt is past directed 
so that J""*" is given by i = and r — > +00, respectively. 

We note that in our formulation of the equations, choice of gauge and with the 
choice of data above, the following "boundary conditions" for the fields Eat and Bab 
are implied at each instant of time 'W . Introducing the new fields 

-^11 '■— (-^11 + ^22)12, E22 := {Ell — £'22)/2, 

and similar for the magnetic part i?ab, we find on both symmetry axes at each instant 
of time 

£'12 = i?13 = E22 — E23 = 0, B12 = -B13 = ^22 = ^23 = 0, (8) 

whereas E*-^ and B^^ are not restricted. Similar properties are obeyed by all symmetric 
Gowdy-invariant covariant 2-tensor fields, for instance Xab- Also, one can derive 
analogous conditions for some of the connection cocfhcients which we do not write 
here. In addition, the constraints of the Bianchi system ^ imply on both axes 

dxE*ii = o^e;, = d^Bl, = d^B;^ = O. (9) 
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2.6. Numerical approach 

The numerical method which we apply in this paper will not be discussed here, see 
[5] for that. In summary, it is a pseudospectral single patch approach; as coordinates 
on S^ we decided for the Euler angle parametrization. The relevant singular terms 
are then determined by the coordinate components of the reference frame {Ya} in 
(|6a|) -(|6c |) . We succeeded in regularizing spectrally the formally singular terms in the 
equations [9\. 

For our investigations here we made the following special choices which have 
not yet been fixed in [9J. First, we decided to use no automatic adaption for the 
spatial resolution. Namely, it turned out to be more convenient to control the spatial 
resolution manually in order to study convergence of the code at very late times. As 
the time integrator we decided for the adaptive fifth order "embedded" Runge Kutta 
scheme [35] ■ Further, we found that it is necessary to enforce the boundary conditions 
([5]) at each time step in order to keep the code stable at very late times of the evolution; 
more details about our "partial enforcement" scheme can be found in [9] . Last but not 
least we mention that we did not use the quad-precision option of the Intel compiler 
[28] , Except for very small values of e, see more comments later, double-precision 
turned out to be sufficient and the constraint growth was more or less under control. 

3. Results and their analysis 

3.1. Setup and main results 

Our strategy to analyze the properties of solutions close to the A-Taub-NUT family is 
as follows. For fixed A-Taub-NUT background spacetimes, i.e. for given values of a and 
Eq, we compute and analyze the solutions corresponding to initial data given by several 
small values of e > 0. In [J we have speculated that A-Taub-NUT Cauchy horizons 
might be stable for small e in the sense that a smooth Cauchy horizon eventually 
forms for all < e < Scrit for some Scrit > 0. The analysis there was done on the basis 
of a numerical code, namely the (2 -t- l)-code in [9], which was not well adapted to the 
demands of the problem. Here we report on results obtained with the new (1-1- l)-code 
introduced in [S] which allows much more reliable investigations. The main reason is 
that with the (1 -I- l)-code no spatial resolution is wasted in one of the symmetry 
directions; recall however that a direct reduction to (1 -I- 1) in our formulation is not 
possible. We still use the same gauge and formulation of the equations as in the 
(2 + l)-code. 

Let us choose the parameters a and Eq such that the A-Taub-NUT solution 
corresponding to £ = has a non-degenerate Cauchy horizon, and let txN be the 
corresponding critical coordinate time in Levi-Civita conformal GauB gauge. Our 
investigations suggest the following results which supersede the previous speculations. 

(i) For each e > 0, there is a time < t^ < txN when the curvature blows up; further 

t^ / tTN for e \ 0. 
(ii) The mechanisms driving the blow up are "qualitatively independent" of e (in the 
sense below) ; hence blow up takes place always "in the same way" in this family. 

(iii) The standard continuity property holds; the smaller e > 0, the longer the 
corresponding solutions stays close to the A-Taub-NUT solution given by e = 0. 

These results suggest that SCC is true in this class of solutions. Below we give more 
details on the dynamical processes which seem to give rise to SCC. 
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Figure 2. Curvature blow up for various values of e 

In all of our numerical investigations below we choose a = 0.7 and Eq = 0, and 
the reason is the following. Our aim is to study perturbations of spacetimes with a 
Cauchy horizon, and this is the case for all A-Taub-NUT solutions with, for instance, 
Eq = 0, and a < •\/3/2 « 0.866 (for A = 3). We expect that qualitatively the results 
do not depend strongly on these parameters if the parameters stay in a regime in 
which the background Cauchy horizon is non-degenerate. This is the case e.g. for all 
a < \/3/2 and Eq — 0. The situation might of course be very different when the 
background is of degenerate type, for instance for a = \/3/2 and Eq — 0, however, we 
do not study that case here. For a — 0.7, we have txN ~ 0.69599284. 

figure [2] shows the numerical plot of a spatial nornu]] of the physical Kretschmann 
scalar versus time for various values of £. In that plot, we have normalized the time 
axis such that the critical time txN of the unperturbed solution is at roughly 20. The 
smallest value of e in our numerical series was 10~^, see comments on this below. 
Note that more details on the numerical properties of the solutions, at least in the 
case e = 10^^, are available in [9]. 



3.2. Heuristic discussions of dynamical effects 

3.2.1. Preparations Now we discuss the numerical results in a more detailed way 
and complement them with heuristic arguments in order to give some explanations 
and extrapolations of our results. A rigorous discussion would be clearly preferable, 
however, is believed to be extremely difhcult since many terms with their individual 
sign structures seem to be significant for the dynamics. One advantage of a 
numerical approach is that it can be expected to provide information in particular on 
intermediate evolution times. This might shed light on for instance the question which 
dynamical processes are responsible to implement SCC. The idea for the following 
heuristic discussion is to consider the linearization of the evolution equation for small 
£ > around the A-Taub-NUT solution given by £ = 0. One can believe that this 
yields a good approximation for small e and small t. Unexpectedly, we find below 
that a discussion exclusively based on the linearized equations already fails at early 
times, i.e. in a small past neighborhood of JT"*". However, by taking into account some 
nonlinear effects supplementing the linearized equations, we are able to obtain some 
insights. It should be noted that we are only interested in qualitative arguments now 

1 This norm is the L^-norm on S'^, however, as mentioned in [9], the volume element of the Euler 
angle coordinates is left out in the integral. 



Investigations of solutions of Einstein's field equations close to X-Taub-NUT 13 

and that none of the simphfying assumptions or approximations considered in the 
following are used for the numerical calculations. We point out again, that our time 
coordinate which we always refer to runs backward in time, i.e. when we speak, e.g. of 
the early phase of the evolution we mean the dynamics close to J7^ , hence very far in 
the future from the physical point of view. After these discussions, we also elaborate 
on the limitations of the current results and hence what remains to be done in future 
work. 

In order to proceed with the linearization we write for the unknowns 

'^a "^a '^ ^"^a i Aab A~ab '^ ^Aabi '^ a c '^ a c '^ '^'^ a c^ ■ ■ ■ i 

where higher order terms in e, which are in principle present for some of these 
unknowns, have already been neglected. Below we express some of the background 
A-Taub-NUT components, marked by o, by means of functions ai{t) and 03 (i) and 
their time derivatives such that 

{ej'){t) =diag(ai(i),ai(i),a3W), 

iXabm=diag{-a[{t)/ai{t),~a[{t)/ai{t),-a'^{t)/a3{t)), 
f^\{t) =ai(i)Va3(i)-2a3(t), t,% = -t^^ = -a.itf/asit); 

all other connection coefficients of the background vanish. Here, a prime denotes the 
derivative with respect to the time t of the A-Taub-NUT solution given by e = 0. The 
function 03 (t) equals the initial data parameter a before at i = 0, further ai(0) = 1. 
The explicit form of the functions ai(t) and 03 (i) can be derived from the relations 
given in section [2T3l however, in order to evaluate the transformation from Brill's time 
coordinate r to our time coordinate t one has to compute an integral which in general 
cannot be represented in closed form. In any case, in order to keep the expressions 
below as short as possible, the other three distinct non-vanishing functions of the 
background L^^ = L22, Eii and B^^ are not written out in terms of ai(i), 03 (i) and 
their derivatives in the following. 

Note that for the full evolution equations with our choice of gauge and frame, 
the 1-2- and 2-3-components of all covariant 2-tensor fields in the set of unknowns are 
identically zero for all times at pi = 0. Also some of the components of the connection 
coefficients vanish there. Hence, we do not need to consider their evolution equations. 
Further we restrict our attention now to the evolution equations of the electric and 
magnetic parts of the rescalcd Weyl tensor (the Bianchi system) since those describe 
the evolution of curvature. Only where necessary, also the other equations are taken 
into account. With F{x) ■= — 2cotx, the linearized evolution equations take the 
following form 



" V ' 

Term IV 

3 



?n(fi^ 



2^n(rA 




(10a) 



Term V 
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Term III 




(12!,) 



TermV 

In order to avoid confusions we point out that studying solutions of these 
equations, which originate in a hnearization around a A-Taub-NUT spacetime, is 
fundamentally different from considering the spin-2 system, i.e. the Bianchi system, 
on the same A-Taub-NUT background. Both lead to a similar, but distinct set 
of equations. In the latter case, which is not of interest here, certain consistency 
conditions [37] must be satisfied, otherwise no solutions need to exist at all. 

3.2.2. Dynamical effects at early times 

Oscillatory phase In this earliest phase in the past of J^^ , the "evolution makes 
a decision on the signs of the quantities" by performing a couple of oscillations 
before going over to a more or less monotonic behavior. Consider the linearized 
evolution equations pOap and (|10b[) first. We know from our choice of initial data 
that both £^22 ^^d B22 are zero at i = 0. If this is true at some other time, the 
equations imply that also their time derivatives are zero, if additionally Xn = X22, 
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Fj^^g ~ —^2^3 (term V) and Bi^ and -E13 are proportional to sinx (term IV); note 
here that (29^ + F{x))smx — 0. This is the case initiaUy. Checking the hnearized 
evolution equations for Xii~X22 and ^1^3 + 1^2^31 which are not printed here, one finds 
that also their time derivatives are zero in that case. Further, due to the linearity 
of the equations and due to the required behavior at the boundaries, Ei^ and -B13 
must be proportional to sinx for all times. In summary, the linearized evolution 
equations imply that for these initial data, E22 and B22 are zero identically for all 
times. Nevertheless, we find numerically that E22 and B22 are non-zero for t > 0. 
In fact, they grow so strongly that they become dominant after a short time. Hence 
a discussion which is exclusively based on the linearized equations with these initial 
data would miss this important character of the solution. This fact is related to the 
results in 4 which however restricts to A = 0. A simple nonlinear mechanism driving 
£'22 and B22 away from zero is the following. If due to nonlinear dynamics there is a 
small higher order component, say, 

i?i3(i,x) = sis' (t) sin x + 5^3^ (i) sin 2x (13) 

with b[1^ < 5^3^ , then term IV in (fTUat yields 

i2d^+F{x))Bu^^4Bghm^X- 
Indeed we find numerically that the early dynamics can be described very well as 

E*22{t,x) = E*2i°\t)sm^X; (14) 

the same holds for the magnetic component. Another conceivable lowest order 
contribution to Eig, and Big, is one which is proportional to cos x, but then according 
to the boundary conditions ([8]) , there must also be a term proportional to cos 3x which 
is of higher order and can be neglected so far. Another possibility consistent with the 
boundary conditions is of the form 1 — cos 2x. However then it turns out that term 
IV in pOap and (jlOSp violates ([5]), i.e. the boundary conditions implied by the Bianchi 
constraints. Hence i\l'S\i is the lowest order explanation for (|14p . In any case, as soon 
as E22 and B22 have been driven away from zero, the linearized equations govern their 
evolution, in particular terms I and II, see below, term IV does not become significant 
before B[^' in (fT51) has become large enough in comparison to BiJ; the same for 
the electric component. At early times, we also find that terms V can be neglected. 
One explanation is the following. When t is small, the only significant sources in the 
linearized evolution equations for the quantities xii ~ X22 and ^1^3 + r2^3 ^^^ ^^e 
quantities E22 and -B22 which are however very small. Hence these quantities are even 
smaller. Now, since E^i and Bl^ are not larger than order unity for all times, this 
implies that term V can be neglected for small t. What is the effect of terms I and II? 
At early times ai and 03 differ only slightly, further a^ah ~ 1 (6 = 1,3) due to the 
choice of initial data, and thus term II can be of the same order of magnitude as 
term I. Indeed, the sign structure of these terms implies oscillations of the amplitudes 
£^22 (i) and ,622 V-) i^ time. Now, both of ai and 03 are monotonically increasing, 
but 03 eventually stronger than ai. These oscillations stop when 03 becomes so big 
in comparison to oi that term II cannot have an significant effect anymore; note that 
in particular a'^/a^ is strongly monotonically increasing as being proportional to xaa- 
The oscillations die out relative early, and in our notation above E22 {t) < and 
-^22 (^) ^ ^ from then on. Then, term I leads to strong growth of the absolute values 
of the amplitudes. We find numerically, that [£'22 (01 > 1^22 (01 for all following 
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early times so that there is a net effect for the Kretschraann scalar Q. Later in the 
evolution, this growth becomes dominant and eventually seems to lead to blow up of 
the Kretschmann scalar. In figure [31 we see a plot of the spatial distribution of these 
components at a time for which our discussion so far is still valid. However, first signs 
of further nonlinear effects become significant which we will discuss later. 

Next consider the linearized evolution equations (|lla|) and (jllb|) of E13 and B13 . 
The early evolution takes place in a similar way as before. The main difference is 
now however that the initial data of E13 is non- vanishing and hence the linear effects 
described by these equations dominate the nonlinear ones at the beginning. In the 
same way as above, terms I and II are dominant with a correction in the oscillatory 
structure by the additional term III. Again, oscillations stop when term II becomes 
insignificant in comparison to term I. The infiuence of term V is different than before. 
Namely, by checking the signs and magnitudes of Fg^g and X13 numerically, term V 
can be seen to become significant a bit later in the evolution. Then it decelerates the 
growth of Ei3 and B13, and so they do not blow up as quickly as £'22 ^^"^ ^22 above. 
The significance of term IV becomes obvious later. For the whole early time, -E13 and 
Bi3 can be described very well as proportional to sinx, and we find numerically that 
their amplitudes are negative from the time on at which oscillations stop; consider 
again figure [31 In particular, the plot shows nicely the different rates of growth for 
the components discussed so far; compare also to figure [H 

The discussion of E^i and Bl^ and their linearized evolution equations (|12ap and 
(|126|l is very similar. Again term I and II are dominant at early times and lead to 
oscillations in time with corrections by term III. Eventually, after these have stopped, 
term V decelerates the growth similarly as in the 1-3-case so that they grow much 
less than in particular E22 and B22- Both E^-^ and B^^ can be described well as 
proportional to cosx with eventually positive amplitudes. 

For understanding the limit e ^ 0, one should be worried whether the relevant 
terms in the equations are of different relative magnitude for different values of e. If 
yes, there would be a dependence of the character of the early oscillations on e, and 
this could lead to a different sign structure of the solution eventually. Most likely 
this would mean that the late time behavior, i.e. close to the critical time, depends 
strongly on e. However, the structure of the linearized evolution equation and the 
initial data are such, that the amplitudes of the dominant terms at early times have 
a consistent scaling with e. Hence the early oscillations are, apart from this scaling, 
independent of e. In detail, the dependence of the initial data on e suggests 

EiitB*i,Ei3,Bi3 ^ e, (15a) 

which is compatible with the equations. Further it is reasonable to expect 

i^2*2,^2*2~e2_ (15^) 

Indeed, one can check numerically that this scaling holds true for early times, but we 
will not show this here. As it can be expected, the behavior at later times cannot be 
described in this simple way anymore. 

Intermediate nonlinear phase We proceed with the discussion of that phase of the 
evolution in which nonlinear effects take over. In that phase, the evolution forms the 
spatial shapes of the unknowns, and leading order effects are the following. As we have 
said before, we find p^ for early times induced by the higher order component in (|T3]); 
similarly for the magnetic component. From the sign structure of £^22 and i?22 '^^ can 
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Figure 3. Spatial depen- 
dence at t = 0.69598373 for 
e = 10-4 



Figure 4. Spatial depen- 
dence at t = 0.69599118 for 
e = 10-* 



read off that the sign of this higher order component must be such that the spatial 
extrema of E13 and B13 slowly move to the south pole. Further we find that term IV 
in (|llap and (jllbp is then also proportional to sin 2% with the "right" signs, hence this 
effect is amplified. In fact this is what we observe numerically, see figure |31 Certainly, 
such nonlinear interactions are also able to produce a component proportional to sin 3x 
or (cosx — cos3x), recall the boundary conditions ([5]), for £^13 and B13. For the latter 
contribution, term IV in (jlOap and (jlO&P yields a violation of (O, hence must not be 
present. Consider sin 3%. With this, term IV in (|10ap and (jlO&P obtains a component 
proportional to cosxsin^ x- The effect is that also the extrema of £'22 and ^22 start 
to move slowly, however, with this simple description we are not able to deduce the 
signs of that term and hence the direction of the move. Numerically we see that 
these extrema move toward the north pole in figure [3] and [H Now, considering again 
term IV in (jllal) and (jllSp with this new component proportional to cos x sin x, the 
quantities £13 and Bi^ obtain an additional contribution of the form sin3x, however, 
with an amplitude opposite to the previous one; thus the nonlinear interactions must 
dominate the linear ones here to keep this process alive. We find numerically that the 
signs eventually become such that the quantities E13 and B13 grow positively close to 
the north pole, while the move of the extrema to the south pole is decelerated. All 
this is visible in figure [4l Indeed we find eventually that the positive growth of £13 
and Bi3 at the north pole becomes much stronger than the dynamics at the south 
pole. 

When comparing different values of e, one finds that the smaller e > 0, the 
later, but eventually the more rapid, all these nonlinear effects take place. One may 
certainly wonder what is special about, say, the north pole in these spacetimes. This 
issue becomes clear when one recalls that going from s to —s is equivalent to the 
interchange of the north and south pole. 



3.2.3. Universality and the late time behavior The longer the evolution proceeds, 
the higher is the order of the significant nonlinear effects. We make no effort here 
to describe all those. Nevertheless, we find numerically that the qualitative picture 
described before does not change anymore for late times. Some of the quantities tend 
to blow up and there is an even stronger concentration of curvature close to x = 0. 
Further, the computations suggest that qualitatively the same picture is true for all 
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small values of e. Certainly, nonlinear effects become significant at a later time the 
smaller e is. In particular the magnitude of E22 and B22 is smaller in relation to the 
other components at early times according to the scaling given by (|15al) and (|15&l) . But 
from some time on, and then on a shorter time scale, these overtake and a "universal" 
qualitative picture is obtained. All this is demonstrated in the numerical plots in 
figure [5] to m As already mentioned, we see particularly strong blow up of E22 and 
B22- Indeed, this seems to yields the main contribution to the Kretschmann scalar 
at late times. According to the plots above, Ei^ and Bi^ also grow at the north pole 
and maybe even blow up eventually. However, this growth seems to be slower. Why 
does the blow up of curvature seem to be confined to a small region around the north 
pole? We can give an explanation for this below, and this is most likely caused by our 
choice of foliation. 

However, what happens exactly on the north pole? The only curvature 
components which do not vanish there according to ([5]) are Eii and B^^ which 
we have said only little about so far. It seems that, as long as we have evolved 
numerically so far, those stay bounded and the dynamics is very slow. According to 
their evolution equations, their behavior at x = depends strongly on the behavior 
of the x-derivatives of E13 and i?i3 there. The numerical results suggest that those 
are positive and grow, but relatively slowly. In fact, it seems that their growth is even 
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slower in comparison to the growth of E22 and B22 the smaller e is, see figure [5] to 
[5] In total, it must be said that it is so far an outstanding question if the curvature 
at X = stays bounded, in particular for small values of e. In order to demonstrate 
how difficult these investigations are with the current numerical setup, figure [9] shows 
the behavior of E^^ at x = for times shortly before the numerical solutions blow 
up for various final numerical resolutions. Here N is the final number of spatial grid 
points in the x-direction and h is the final value of the size of the time step, i.e. 
equals hmin as given in [0]. In addition, figure [TU] shows the behavior of the constraint 

violations Norm'^'^"""' '^'' defined in [9]. This plot suggests that the numerical solution 
converges and that the constraint error can be kept under control for sufRciently high, 
in particular spatial resolutions at least temporarily. Nevertheless, since on the one 
hand the demand for resolution at those late times becomes so high, but on the other 
hand the dynamics of E^^ (also B^^) at x = seems to be so slow, we cannot answer 
the question about the behavior of curvature on the north pole yet. This means that 
we are also not able yet to exclude the possibility that a strong blow up of El^ and 
B^i at very late times induces a drastic change of the whole asymptotic dynamics. In 
any case, ideas how to optimize the numerical late time behavior and the setup of our 
approach are discussed below. 



3.3. Further comments on the solutions 

3.3.1. Area of the group orbits In section [2.41 we have introduced the notion of the 
rescaled orbit area for the S'^-Gowdy class and have stated that one can believe that 
a relevant subset of the MGHDs can be covered by areal foliations. In the T'^-Gowdy 
case without cosmological constant the singularity is attained generically exactly where 
the group orbit area density vanishes. If we are willing to assume that this is true 
also for our S'^-Gowdy solutions with A > 0, then we can get a feeling on how "close" 
the numerical evolutions have eventually approached the "singularity" at the end 
of the numerical runs by looking at the rescaled orbit area density. Gertainly, our 
foliation is not of areal type, hence this quantity is not constant on our t = const 
slices, see figure [TT] for a typical late time profile. According to this interpretation, 
we see that our foliation has approached the singularity further at the north pole 
than at the south pole which explains why curvature blows up much more strongly 
there. In fact, figure [T2l implies that this discrepancy between the poles, where the 
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Figure 12. Extremal values 
of rescaled orbit area for e = 
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on t = const-slices 



maximum and minimum respectively of this rescaled quantity are always attained in 
our case, becomes larger the further the evolution proceeds. Thus the singularity is 
approached "inhomogeneously" in this sense in our gauge. We just note that this is 
not surprising since many gauges which have been used in the literature so far have 
such a property. Hence one can expect that only special gauges are suitable for our 
kinds of investigations, and it will be a matter of future research to experiment with 
such other foliations, see below. Concerning the question of SCO our gauge is thus not 
appropriate enough. For other open questions, however, as that about the appearance 
of Gowdy spikes, our computations have the potential to yield new insights. 



3.3.2. Gowdy spikes From the results obtained in the T'^-Gowdy case, one can expect 
spikes in the sense of section 12.41 to form under certain conditions at isolated points 
away from the axes. The previous discussion about the gauge however suggests that 
we cannot conclude that the local curvature blow up close to the north pole is such a 
spike. Rather, it can be caused by the "inhomogeneous" approach to the singularity. 
It is nevertheless possible that spikes show up later in the evolution even in our gauge, 
but most likely they would be difhcult to distinguish from our gauge effects. However, 
from our numerical data obtained in Levi-Civita conformal Gaufi gauge we can - in 
principle - compute the geometrical quantities with respect to an areal foliation until 
very "close to the singularity" in a small spatial neighborhood of the north pole. With 
this it would be feasible to study the possible appearance of spikes there. However, 
this has not been done yet. In areal gauge, one often introduces a time coordinate r 
such that the rescaled orbit area is proportional to e^"^, at least for late times in the 
§^-case. Hence, a local reconstruction of the areal foliation from the current numerical 
data would reach Gowdy times of the order r « 20, which is in fact quite impressive. 
Certainly, it is possible that our family of solutions given by low e corresponds to 
"low-velocity" solutions. The concept of the asymptotic velocity is not defined here, 
see [29l |45] for the S^-case. In fact the relation of this asymptotic quantity at the 
singularity and our definition of £ as a data parameter on >/+ is not easy to see. It 
sounds reasonable that the larger e is, the higher the asymptotic value of the velocity 
becomes. Now, in the T^-Gowdy case in vacuum with A = 0, one finds that spikes 
occur only if the asymptotic velocity exceeds a certain value. So we would expect 
spikes to form in our class of spacetimes away from the axes only for large e. It is 
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an unexpected outcome of our numerical experiments that it seems easier to compute 
solutions for larger e than for small values, see below. Hence there is hope that spikes 
can be studied with the current numerical setup. This is planned as a future research 
topic, but in this paper we are rather concerned with the limit e —f 0. Nevertheless, 
we would like to note that our numerical method can be expected to be able to resolve 
Gowdy spikes in principle because of our experience with numerical investigation of 
T'^-Gowdy solutions in ^. The question of the behavior of the solutions on the axes 
is particularly interesting since, as discussed in [45], although one can expect some 
discontinuous behavior at the axes in a particular family of parametrizations of the 
metric, it is not clear whether generic S'^-Gowdy solutions have a true spike there. 
Thus it is remarkable that our numerical results suggest that even at those relatively 
late Gowdy times, the curvature has not blown up on the axes yet. However, this 
might of course also be a consequence of the special choice of initial data. 

3.4- Limitations and prospects for future work 

Our approach of studying the properties of this class of solutions has several limitations 
so far, some of which were already mentioned above. Some of these limitations are 
caused by our non-rigorous combination of numerical and heuristic analyzes, and are 
of principle nature. Others, however, can be avoided more easily and are currently 
work in progress. In fact we believe that the numerical method has not been pushed 
to its limits yet. 

In our investigations, we are not yet able to conclude finally on the question 
of sec in this class of spacetimes, although the numerical results suggest so far 
consistently that SCC holds in the sense of section 13.11 Nevertheless, we cannot 
yet deduce the very late time behavior of curvature at the north pole, mainly, but 
not exclusively, due to numerical problems. Also away from the north pole, we are 
facing the problem that the foliation does not proceed toward the "singularity" fast 
enough as was discussed above. Concerning this gauge issue, we would hence like 
to try various other gauge conditions in future work, for instance areal gauge, which 
might however be problematic for A > in general, but possibly not here. Within 
the class of general conformal Gaufi gauges there is also freedom left which we have 
not exploited yet. Nevertheless, it is expected that the description of the transport 
of spatial symmetries is not trivial in general in this gauge. We would like to note 
that there is also gauge freedom left within the Levi-Civita conformal GauB gauge, 
namely to perform an arbitrary conformal transformation of the initial 3-metric and 
correspondingly of the other data components. It has not been investigated yet if 
it would be possible to optimize the foliation in this simple way. Other promising 
gauges are those of harmonic type and constant mean curvature foliations since there 
is some experience with these in the numerical relativity community. However, for 
their use with the conformal field equations, it must be investigated first if and how 
those gauges can be put into the required form. 

Concerning the numerical method, there are currently a number of obstacles for 
driving the numerical solutions "further to the singularity" . Some of those are straight 
forward to solve, others are not. The runs presented here each took a few hours on 
single processor machines. It is of course no principle problem to let the runs continue 
for a few days etc. But at the moment, we have a problem with disk space since our 
data output routines are not yet optimized; clearly, this problem is straight forward to 
fix. What are the main reasons that it is numerically demanding to continue the runs? 
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At the north pole, the solutions develop large gradients, hence much spatial resolution 
is needed there. So far in our code, resolution cannot be increased locally, and thus 
increasing the spatial resolution to cope with the demands at the north pole always 
goes hand in hand with a waste of grid points around the south pole where relatively 
little happens. Up to now, we have investigated only little on spatial coordinate 
transformations which would put more spatial points close to the left axis and less to 
the right axis without a direct change of the numerical infrastructure. In any case, 
this trick will be studied further in the future. It should be noted that our code is 
not yet optimized technically for high spatial resolutions in general; in particular we 
still do not use the Fast Fourier Transform (FFT) [15] yet but only partial summation 
[10] . Also it may be true that there are more optimal time integrator schemes than the 
adaptive Runge Kutta schemes of our choice, comments on this can be found in [TO]. 
Further, it might make sense to think about parallclization of the code; this should 
be straight forward with some publically available FFT libraries as for instance |23j . 
It is certainly a justified question if a pseudospectral method like ours is suitable at 
all for this class of problems. Namely, although these methods are highly accurate 
for lower resolutions they might be too slow for high resolutions. Thus it makes 
sense to also investigate into other methods, for instance finite differencing methods, 
maybe even with multipatch or mesh refinement. For instance the similar numerical 
investigations of Gowdy spacetimes with spatial §^ x §^-topology in [23] were done 
with finite differencing. Another non-related numerical issue shows up for small values 
of e at late times. Then, the numerical noise caused by round-off errors can become 
significant. Indeed, this was the main reason why we did not decrease e further than 
10~^. Maybe, the only solution to this problem is to switch to quad-precision which is 
possible for Intel compilers and Intel processors, but then the code runs much slower. 
Further comments on our numerical infrastructure can be found in [9]. 

A further aspect one should keep in mind is the question about the optimal 
formulation of Einstein's field equations for our purposes. It is justified to ask if it 
is necessary to use the conformal field equations for our applications at all. There 
are two motivations why we decided for these equations. First, prescribing data on 
J'~^ is simpler than in the case of Cauchy surfaces of the physical spacetime. Second, 
we wish to study FAdS spacetimes since those are motivated by current cosmological 
observations and, which is maybe more important here, the MGHD of any data on 
J'~^ is geodesically complete to the future by construction. So SCC is obeyed at 
least to the future. In contrast, prescribing data on a Cauchy surface of the physical 
spacetime yields no a priori knowledge in general whether the corresponding solution 
is future geodesically complete. Of course, one might want to restrict to the question 
of SCC in the past direction as a first step. Then, the conformal field equations can 
be avoided. It is possibly an advantage to rely on the experience of the numerical 
relativity community which is available for certain formulations of EFE but not so 
much for the conformal field equations. For instance there exist methods to damp 
constraint violations etc. Indeed, we have performed successful numerical tests with a 
well-adapted system other than the CFE in the T^-Gowdy case in [7] and would like 
to work out a similar system for the S^-case. 

4. Summary 

In this work, we have presented investigations of the evolution properties of a special 
class of spacetimes given by certain perturbations of the A-Taub-NUT solutions. By 
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means of a combination of numerical and heuristic analyzes, we were able to draw 
a qualitative picture of the late time behavior of the solutions. The results suggest 
that sec holds in this family and we are able to give more details on the dynamical 
processes which seem to give rise to SCC. However, we are not yet able to make final 
conclusions. Nevertheless, we identify some of the problems and limitations of our 
current approach and give prospects for future work. It is surprising that the issue of 
strong cosmic censorship and of the BKL-conjecture is still open even in this relatively 
simple class of §^-Gowdy spacetimes. Certainly the family of spacetimes considered 
here is very special, and it is unclear what information we yield both about the general 
SCC conjecture and about SCC in the class of general S'^-Gowdy spacetimes. We 
believe however that the results presented here are promising steps in order to derive 
a consistent picture on these questions in the near future. 
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